Marine heatwaves alter the nursery function of coastal habitats for juvenile Gulf of Alaska Pacific cod

Marine Heatwaves (MHWs) can directly influence survival of marine fishes, particularly for early life stages, including age-0 juveniles during their residence in coastal nursery habitats. However, the ability of nurseries to support high fish densities, optimize foraging and growth, and protect against predators may be altered during MHWs. Gulf of Alaska Pacific cod (Gadus macrocephalus) larval, juvenile, and adult abundances declined dramatically following MHW events in 2014–2016 and 2019. To evaluate coastal nursery function during MHWs, we compared diet composition, recent growth, size, condition, and abundance of age-0 juveniles throughout their first summer before, during, and between MHWs. Diet shifted to larger prey during MHWs, particularly mysids, but diet did not appear to influence growth. We observed faster growth rates during MHWs, yet even when accounting for growth, we could not explain the higher body sizes observed in August during MHWs. Together with lower abundance and the near absence of small fish in the nursery by August during MHWs, these patterns highlight potential for size-selection and a reduced ability of nursery habitats to buffer against environmental variability during MHWs, with only a small number of large “super survivors” persisting through the summer.

Table S2: Sample sizes for total catch, diet composition, and growth analyses of juvenile Pacific cod by year for July and August sampling.Summary statistics are reported for Trident Bay, AK temperature (ºC), age-0 Pacific cod catch per unit effort (CPUE), standard length (mm), and body mass (g).Samples were not available for July 2006July , 2008July , and 2011, or , or               Table S6: Steps for growth model selection, using a three-step process: 1: Determine optimal random structure; 2: Determine optimal fixed structure; 3: Determine optimal error structure.All predictor variables were scaled and centered.

July August
Step  Step 3: Optimal Error Structure (calculated with REML) Step

Figure S1 :
Figure S1: Water temperatures from Trident Bay, Kodiak, AK were >3C warmer during marine heatwave events than temperatures prior to the heatwaves.Water temperatures between the heatwaves (2017-2018) remained higher than normal, although not as warm as the heatwaves.Temperature throughout the year is shown by heatwave class (left) and by year (right) from 2006 to 2019.Shaded polygons represent 95% confidence intervals.Water samples are collected from 10m below mean lower low water.

Figure S2 :
Figure S2: Juvenile Pacific cod body mass (a) and length-weight condition residuals (b) in July (light grey) and August (dark grey) from 2006-2019.Plots are shaded based on their heatwave class, with blue representing years before the MHWs, red representing years during MHWs, and light grey representing years between MHWs.Boxplots extend from the first to the third quartiles of the data, with whiskers that extend to the largest values no further than the 1.5*IQR.Individually plotted points represent outliers beyond the 1.5*IQR range.

Figure S3 :
Figure S3: Juvenile Pacific cod square root-transformed stomach fullness in July and August did not vary across marine heatwave events.Light-grey boxplots represent July stomach fullness, and dark-grey boxplots represent August stomach fullness.Plots are shaded based on their heatwave class, with blue representing years before the MHWs, red representing years during MHWs, and grey representing years between MHWs.Boxplots extend from the first to the third quartiles of the data, with whiskers that extend to the largest values no further than the 1.5*IQR.Individually plotted points represent outliers beyond the 1.5*IQR range.

Figure S4 :
Figure S4: Juvenile Pacific cod diet composition (represented by year) varied by year and by heatwave class, with mysids representing a large portion of the diet by August in marine heatwave years.Diet composition is represented by percent Prey Specific Index of Relative Importance (PSIRI) between 2006 and 2019 for July (top panel) and August (lower panel).The "other prey" category includes all prey that was not found in at least 3.5% of all stomachs sampled.Years are colored by their heatwave class, with blue representing years before marine heatwaves, red representing years during marine heatwaves, and grey representing years between marine heatwaves.

Figure S5 :
Figure S5: Size-restricted nonmetric multidimensional scaling (NMS) biplot of July (right) and August (left) diet composition for juvenile Pacific cod from 2007-2019 in heatwave conditions (red), before heatwaves (blue), and between heatwaves (grey).Sizes were restricted to 50-60 mm standard length in July and 70-90 mm standard length in August.Ellipses represent 95% confidence intervals in the two major axes.

Figure S6 :
Figure S6:Predicted August size distributions (beige) and actual August size distributions (grey) for (a) real data collected from the field; (b) six simulated populations of 10,000 individuals (one each for predicted August sizes and observed August sizes across the three MHW classes); and (c) random draws of smaller subsets of the large, simulated population.Predicted August sizes were based on July sizes and July growth rates and calculated by summing mean July growth rates for that year over the number of days between July and August sampling in that year.

Figure S7 :
Figure S7: Juvenile Pacific cod were sampled via beach seine surveys from Anton Larsen Bay and Cook Bay on the northeastern coast of Kodiak Island, Alaska, U.S.A. in July and August 2006-2019.Daily nearshore temperature readings were taken from long-term temperature loggers in Trident Bay on the northeastern coast of Kodiak Island.

Figure S8 :
Figure S8: Relationship between juvenile Pacific cod log-transformed whole body wet weight (g) and log-transformed standard length (mm) from both July and August between 2006 and 2019.

Figure S9 :
Figure S9: Top panel: scree plots for July (a) and August (b) assessing optimal number of dimensions for Non-Metric Multidimensional Scaling ordination based on overall stress.Middle panel: Shepard plots for July (c) and August (d) evaluating ordination fit.The observed dissimilarities within the data are plotted against the Bray-Curtis distances of the ordination.Bottom panel: goodness of fit for July (e) and August (f) to determine the proportion of overall variance explained by each point.Larger circles represent individual Pacific cod stomachs that account for a larger proportion of the overall variance within the ordination.

Figure S11 :
Figure S11: Mean daily relative growth trajectories (mm/mm/day); back-calculated from daily otolith increment widths using a Biological Intercept for the final three weeks of nearshore nursery residence by year (top) and heatwave class (bottom) for July (left) and August (right).Colored polygons represent 95% confidence intervals.Increment number 0 represents growth 21 days prior to capture.Increment number 21 represents growth on the final day prior to capture.Only years with data for both months are shown here.

Step 2 :
Optimal Fixed Structure (calculated with Maximum Likelihood) log(Growth) ~ Size + Temperature + Heatwave Class -2027.85log(Growth) ~ Size x Temperature x Heatwave Class -2550.42log(Growth) ~ Size x Temperature x Heatwave Class -2023.63 log(Growth) ~ Size x Heatwave Class -2526.13log(Growth) ~ Size x Temperature -2022.67 log(Growth) ~ Size x Temperature -2515.12log(Growth) ~ Temperature x Heatwave Class -1990.72 log(Growth) ~ Size + Temperature + Heatwave Class -) ~ Size x Temperature x Stomach Fullness x Heatwave Class -1291.20 log(Growth) ~ Size x Temperature x Stomach Fullness x Heatwave Class -1662.79log(Growth) ~ Size x Temperature x NMS Axis 1 Score x Stomach Fullness x Heatwave Class -1272.07log(Growth) ~ Size x Temperature x NMS Axis 1 Score x Stomach Fullness x Heatwave Class -1612.71log(Growth) ~ Size x Temperature x NMS Axis 1 Score x Stomach Fullness x Abundance x Heatwave Class -1264.71log(Growth) ~ Size x Temperature x NMS Axis 1 Score x Stomach Fullness x Abundance x Heatwave Class -1588.65 August 2011 and 2016.Summary statistics are reported as mean ± standard error.

Table S4 :
Type III SS ANOVA results for evaluating the impact of heatwave (before, during, between) and grouping variable (actual August size and predicted August size) on August standard length.

Table S5 :
July and August prey lists.

Table S7 :
Marginal mean effects for July and August relative growth (mm/mm/day) for best-fit models before, between, and during marine heatwaves.All values are adjusted to otolith increment number 11.This is the companion table for Figure5in the main text.